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13.  ABSTRACT  '  Maximum  JCO  words) 

Robotic  manipulators  will  play  a  significant  role  in  the  maintenance  and  repair  of  space  stations  and  satellites,  and  other 
future  space  missions.  Robot  path  planning  and  control  for  the  above  applications  should  be  optimum,  since  am 
inefficiency  in  the  planning  may  considerably  nsk  the  success  of  the  space  mission.  This  paper  presents  a  global 
optimum  path  planning  scheme  for  redundant  space  robotic  manipulators  to  be  used  in  such  missions.  In  this 
formulation,  a  variational  approach  is  used  to  minimize  the  objective  functional.  Two  optimum  path  planning  problems 
are  considered:  first,  given  the  end-effector  trajectory,  find  the  optimum  trajectories  of  the  joints,  and  second,  given  the 
terminal  conditions  of  the  end-effector,  find  the  optimum  trajectories  for  the  end-effector  and  the  joints.  It  is  explicitly 
assumed  that  the  gravity  is  zero  in,  and  the  robotic  manipulator  is  mounted  on  a  completely  free-flying  base  (spacecraft) 
and  the  attitude  control  (reaction  wheels  or  thrust  jets)  is  off.  Linear  and  angular  momentum  conditions  for  this  system 
lead  to  a  set  of  mixed  holonomic  an  nonholonomic  constraints.  These  equations  are  adjoined  to  the  objective  functional 
using  a  Lagrange  multiplier  technique.  The  formulation  leads  to  a  system  of  Differential  and  Algebraic  Equations 
(DAEs)  and  a  set  of  terminal  conditions.  A  numerical  scheme  is  presented  for  forward  integration  of  the  above  system  of 
DAEs,  and  an  iterative  shooting  method  is  used  to  satisfy  the  terminal  conditions.  This  approach  is  significant  since 
most  space  robots  that  have  been  developed  so  far  are  redundant.  The  kinematic  redundancy  of  space  robots  offers 
efficient  control  and  provides  the  necessary  dexterity  for  extra-vehicular  activity  or  avoidance  of  potential  obstacles  in 
space  stations. 
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ABSTRACT 


Robotic  manipulators  will  play  a  significant  role  in  the  maintenance  and  repair 
of  space  stations  and  satellites,  and  other  future  space  missions.  Robot  path  plan¬ 
ning  and  control  for  the  above  applications  should  be  optimum,  since  any  inefficiency 
in  the  planning  may  considerably  risk  the  success  of  the  space  mission.  This  paper 
presents  a  global  optimum  path  planning  scheme  for  redundant  space  robotic  ma¬ 
nipulators  to  be  used  in  such  missions.  In  this  formulation,  a  variational  approach 
is  used  to  minimize  the  objective  functional.  Two  optimum  path  planning  problems 
are  considered:  first,  gi’  en  the  end-effector  trajectory,  find  the  optimum  trajectories 
of  the  joints,  and  second,  given  the  terminal  conditions  of  the  end-effector,  find  the 
optimum  trajectories  for  the  end-effector  and  the  joints.  It  is  explicitly  assumed  that 
the  gravity  is  zero  ,  and  the  robotic  manipulator  is  mounted  on  a  completely  free- 
flying  base  (spacecraft)  and  the  attitude  control  (reaction  wheels  or  thrust  jets)  is 
off.  Linear  and  angular  momentum  conditions  for  this  system  lead  to  a  set  of  mixed 
holonomic  and  nonholonomic  constraints.  These  equations  are  adjoined  to  the  ob¬ 
jective  functional  using  a  Lagrange  multiplier  technique.  The  formulation  leads  to  a 
system  of  Differential  and  Algebraic  Equations  (DAEs)  and  a  set  of  terminal  condi¬ 
tions.  A  numerical  scheme  is  presented  for  forward  integration  of  the  above  system 
of  DAEs,  and  an  iterative  shooting  method  is  used  to  satisfy  the  terminal  conditions. 
This  approach  is  significant  since  most  space  robots  that  have  been  developed  so  far 
are  redundant.  The  kinematic  redundancy  of  space  robots  offers  efficient  control  and 
provides  the  necessary  dexterity  for  extra-vehicular  activity  or  avoidance  of  potential 
obstacles  in  space  stations. 
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1.  Introduction 


Space  exploration  is  a  new  frontier  in  current  science  and  engineering  [1].  Benefits 
from  the  space  explorat:on  are  enormous,  however,  the  stake  is  also  very  high.  Space 
missions  are  hazardous  to  astronauts  [2]  because  of  extremes  of  temperature  and 
glare,  and  possible  high  level  of  radiation.  The  extra- vehicular  activity  also  consumes 
considerable  time  and  may  need  the  dexterity  and  high  load  handling  capacity  that 
astronauts  can  not  prov  de.  Therefore,  using  robots  in  space  is  beneficial  for  extending 
on-orbit  time  of  space  shuttle  and  increasing  productivity  of  space  mission. 

The  use  of  robotic  manipulators  in  space  applications  introduces  several  new  prob¬ 
lems  which  do  not  arise  in  ground  base  robot  applications.  For  example,  the  motion 
of  a  space  manipulator  can  cause  the  base  (satellite  or  spacecraft)  of  the  manipulator 
to  move  and  disturb  the  trajectory  of  the  spacecraft  [1-4].  This  can  severely  affect  the 
spacecraft  performance  specially  when  the  mass  and  the  moment  of  inertia  of  the  ma¬ 
nipulator  arms  and  the  payloads  are  not  negligible  in  comparison  to  the  manipulator 
base. 

One  solution  to  the  above  problem  is  to  use  reaction  (thrust)  jets  to  control  the 
attitude  of  the  spacecraft  (or  the  robot  base)  [5-6].  To  meet  this  requirement,  the 
spacecraft  must  carry  additional  reaction  jet  fuel.  Since  a  spacecraft  can  carry  only 
a  limited  load,  this  aporoach  may  force  removal  of  other  facilities  of  considerable 
importance  for  the  success  of  the  mission.  Thus  it  may  alter  (reduce)  the  goal  of 
the  mission.  Furthermore,  exhaust  from  the  reaction  jets  may  interfere  with  proper 
operation  of  the  instruments  on  the  board.  For  example,  the  exhaust  may  reduce 
the  vision  distance  of  a  camera  if  the  camera  must  see  along  the  reaction  jets,  and  in 
some  extreme  cases  the  deposition  of  the  exhaust  on  the  camera  lenses  may  completely 
block  the  vision  of  the  camera. 

Another  solution  to  the  above  problem  is  to  keep  the  reaction  jets  off  and  to 
move  the  space  manipulator  arms  in  such  a  fashion  that  it  accomplishes  the  desired 
space  tasks  and  yet  maintains  the  stability  and  the  overall  trajectory  of  the  space¬ 
craft  and/or  the  satellite.  This  approach  will  considerably  reduce  the  reaction  jet 
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fuel  needed  for  attitude  control  and  increase  the  life  span  of  the  spacecraft  and  the 
manipulator. 

It  is  clear  that  the  second  approach  is  far  superior  than  the  first.  For  this  reason, 
most  investigators  in  this  area  have  focussed  their  research  interests  on  the  second 
approach  [1-10].  Lindberg,  Longman,  and  Zedd  [1]  address  various  issues  related  to 
free-flying  space  manipulators  and  provide  a  comprehensive  review  of  several  papers 
on  the  subject.  Umetan  and  Yoshida  [7]  present  a  Jacobian  matrix  formulation  for  the 
study  of  kinematics  and  control  of  a  free-flying  space  manipulator.  Their  formulation 
includes  both  the  linear  and  the  angular  momentum  conservation  conditions.  It 
should  be  noted  that  the  pseudo  Jacobian  inverse  solution  of  the  velocity  constraint 
equations  provides  only  local  minimum  of  the  generalized  velocity  norm  which  may 
not  be  the  global  minimum. 

Recently,  Vafa  and  Dubowsky  [2]  have  presented  the  concept  of  a  Virtual  Manipu¬ 
lator  (VM)  for  space  robots  and  have  shown  that  a  space  robot  and  the  corresponding 
virtual  manipulator  give  identical  kinematic  response.  They  have  also  shown  that  if 
the  net  linear  momentum  of  the  space  robotic  system  is  zero,  then  the  location  of 
the  virtual  ground  remains  fixed.  Thus,  using  the  virtual  manipulator  concept,  the 
studies  and  results  for  ground  based  manipulators  may  be  extended  to  space  manip¬ 
ulator  systems.  Papadcpoulos  and  Dubowsky  [8]  have  used  this  concept  to  study  the 
singularity  of  space  manipulators. 

Longman,  Lindberg  and  Zedd  [3]  present  a  reaction  wheel  compensation  method 
for  attitude  control  of  the  spacecraft  (or  the  space  manipulator  base).  Since  the 
reaction  wheels  use  photo-voltaic  energy,  this  method  uses  considerably  less  control 
jet  fuel  than  a  reaction  jet  control  method.  Reaction  wheels,  however,  do  not  control 
the  translational  disturbances  of  the  spacecraft.  Reference  [9]  provides  the  kinematics 
and  workspace  analyses  of  a  satellite-mounted  robotic  system. 

References  [1-9]  largely  discuss  techniques  that  accomplish  certain  tasks  without 
disturbing  the  stability  of  the  systems.  Fernandes,  Gurvits,  and  Li  [10]  present  a 
method  for  near  optimum  attitude  control  of  space  manipulator  using  internal  mo¬ 
tion.  This  formulation  considers  a  two  point  boundary  value  problem  but  it  does 
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not  consider  the  problem  of  the  end-effector  following  a  path.  As  a  result,  the  holo- 
nomic  and  nonholonomic  momentum  conservation  constraints  are  homogeneous  and 
the  solution  for  the  velocity  equations  lie  in  the  null  space  of  the  Jacobian  matrix. 
This  may  not  be  true  if  the  end-effector  must  follow  certain  trajectory.  For  further 
discussions,  readers  may  wish  to  examine  references  [1-10]  and  the  references  therein. 

In  this  paper,  we  present  a  global  optimum  path  planning  scheme  for  a  redundant 
space  robotic  manipulator  flying-freely  in  a  zero  gravity  space  with  reaction  jets  oif. 
Global  here  implies  that  the  formulation  accounts  for  minimization  of  the  functional 
for  the  entire  path  and  rot  for  the  local  point.  The  extra  degrees  of  freedom  of  the  re¬ 
dundant  manipulator  enlarge  the  workspace  of  the  manipulator.  Furthermore,  these 
extra  degrees  of  freedom  may  be  used  to  optimize  certain  functionals  and  avoid  sin¬ 
gularities.  The  formulation  accounts  for  the  holonomic  and  nonholonomic  constraint 
conditions  arising  from  the  momentum  conservation  conditions.  We  consider  the  fol¬ 
lowing  two  optimum  path  planning  problems:  first,  given  the  end-effector  trajectory, 
find  the  optimum  trajectories  of  the  joints,  and  second,  given  the  terminal  conditions 
of  the  end-effector,  find  the  optimum  trajectories  for  the  end-effector  and  the  joints. 
The  formulation  leads  ,o  a  set  of  differential  and  algebraic  equations.  A  numerical 
scheme  to  solve  these  set  of  equations  is  also  presented. 

The  outline  of  the  paper  is  as  follows:  In  section  two,  we  develop  the  momen¬ 
tum  conservation  and  *he  kinematic  conditions  for  the  system.  This  leads  to  a  set 
of  holonomic  and  nonholonomic  constraints.  In  this  section,  we  also  make  some  re¬ 
marks  on  local  optimization  and  show  that  in  this  formulation  the  second  problem 
may  be  considered  as  a  subset  of  the  first  formulation.  In  section  three,  we  use  a 
variational  approach  tc  optimize  an  objective  functional  subjected  to  the  momen¬ 
tum  and  the  kinematic  constraints  developed  in  section  2.  Boundary  conditions  are 
discussed  in  section  four.  In  section  five,  we  present  a  brief  outline  of  the  iterative 
shooting  method.  Section  six  presents  a  numerical  scheme  to  integrate  the  system  of 
differential  and  algebraic  equations.  In  section  seven,  we  consider  a  numerical  exam¬ 
ple  to  show  the  feasibi  ity  of  the  formulation.  Conclusions  are  presented  in  section 
six.  Finally,  in  the  Appendix  we  demonstrate  the  various  matrices  that  appear  in  the 
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formulation  using  a  twc  dimensional  example. 


2.  Momentum  Conservation  and  Kinematic  Conditions 

In  the  mathematical  formulation  for  global  optimum  path  planning  for  space 
robotic  manipulators  to  follow,  we  consider  the  following  two  problems: 

1.  Given  the  trajectory  of  the  end  effector,  find  the  trajectories  of  the  joints  that 
minimize  a  given  objective  functional,  and 

2.  Given  the  terminal  conditions  of  the  end  effector,  find  the  trajectories  of  the 
end  effector  and  the  joints  that  minimize  a  given  objective  functional. 

We  assume  that  the  gravity  is  zero  in  space  and  the  thrust  jets  are  off.  Thus,  there  is 
no  external  force  or  torque  acting  on  the  system  and  the  system  is  freely  floating  in 
the  space.  This  implies  that  the  linear  and  angular  momentum  of  the  system  must 
be  conserved. 

To  develop  momentum  conservation  and  other  kinematic  conditions,  consider  a 
typical  space  robotic  manipulator  consisting  of  na  arms,  one  base,  and  nw  inertia 
wheels  (Figure  1).  As  shown  in  Figure  1,  the  manipulator  base  is  numbered  0,  the 
predecessor  arms  are  numbered  as  1,  2,  •••,  n„  in  increasing  order  and  the  inertia 
wheels  are  numbered  greater  than  na.  The  joint  between  arm  i  and  its  predecessor 
arm  or  the  base  is  called  joint  i.  It  is  assumed  that  each  joint  has  only  one  relative 
Degree  Of  Freedom  (DOF).  Thus,  the  degrees  of  freedom,  nj,  of  a  spatial  and  planner 
space  manipulator  are,  respectively,  (6  +  na  +  nw)  and  (3  +  na  +  nw).  Formulations 
for  a  system  having  joints  with  relative  degrees  of  freedom  more  than  one  may  be 
obtained  in  a  similar  fashion. 

For  simplicity  in  the  discussion  to  follow,  the  manipulator  base,  the  manipulator 
arms,  and  the  reaction  wheels  are  also  called  bodies.  Thus  body  0  represents  the 
base,  bodies  1  to  n„  represent  the  arms,  and  bodies  (na  +  1)  to  ( na  4-  nw)  represent 
the  inertia  wheels.  As  shown  in  Figure  1,  let  r,  be  the  position  vector  of  center  of 
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Arm  na 


Figure  1 .  Schematic  ot  a  Space  Manipulator 
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mass  of  body  i  with  respect  to  (w.r.t.)  an  inertial  frame,  then  the  expressions  for  the 
linear  and  the  angular  momentum  (p  and  L)  of  the  system  may  be  written  as  [11] 

£  m«i‘.  =  P  (1) 

1=0 

Ha 

Y  [lm  +  r{  X  m.f,]  =  L  (2) 

1=0 

where  m,  and  u>i,  are,  respectively,  the  mass  and  the  angular  velocity  of  body  i,  /,  is 
the  moment  of  inertia  of  body  i  about  its  center  of  mass,  and  the  period  on  a  variable 
(.)  denotes  the  total  derivative  of  the  variable  (.)  with  respect  to  time  t. 

In  many  applications,  it  is  more  convenient  to  include  the  mass  and  the  moment 
of  inertia  of  the  inertia  wheels  to  that  of  the  manipulator  base.  Furthermore,  for 
simplicity  in  the  discussion  to  follow,  it  is  assumed  that  the  linear  and  the  angular 
momentum  of  the  system  are  zero,  i.e.  (p  =  L  =  0).  This  is  not  the  limitation  of 
the  formulation.  Formulations  for  p  and  L  not  equal  to  zero  may  be  developed  on  a 
similar  line.  However,  in  this  case  the  formulation  will  be  more  involved.  Following 
the  above  discussion,  equations  (1)  and  (2)  may  be  written  as 


na 

yi  +  m0r0  =  0  (3) 

i=i 

YVm  +  r,  x  m,rf]  +  [I^u 0  +  r'0  x  m^]  +  [I^f]  =  0  (4) 

>=1  J=na+ 1 

where  m'Q  and  /q  are,  respectively,  mass  and  moment  of  inertia  of  the  manipulator 
base  including  the  inertia  properties  of  the  inertia  wheels,  r'0  is  the  new  center  of  mass 
of  the  manipulator  base,  and  u>™b,  ( j  =  na  +  1,  •  •  • ,  na  +  nw)  are  the  angular  velocity 
of  the  reaction  wheels  with  respect  to  the  manipulator  base. 

Equation  (4)  shows  that  the  attitude  of  the  manipulator  base  may  be  controlled 
in  a  relatively  simple  manner  by  controlling  a  set  of  three  orthogonally  placed  inertia 
wheels.  Controlling  the  attitude  of  the  base  without  using  inertia  wheels,  however,  is 
more  difficult. 
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For  a  general  system  one  can  write  a  close  form  solution  for  the  time  integral  of 
equation  (3)  but  not  ol  equation  (4).  Thus,  equations  (3)  and  (4)  provide  a  set  of 
holonomic  (Eq.  3)  and  nonholonomic  (Eq.  4)  conditions.  The  system  must  maintain 
these  momentum  conservation  conditions  which  cause  disturbances  in  the  system. 
The  formulation  to  follcvr  will  show  that  by  proper  momentum  management  one  can 
perform  the  desired  task  and  still  satisfy  the  equations. 

Integrating  equation  (3)  with  respect  to  time,  we  obtain 

53  m<ri  +  moro  =  mrG  (5) 

1=1 

where  m  is  the  total  mass  and  is  the  center  of  mass  of  the  system.  Since  the 
time  integral  of  equation  (3)  must  be  a  constant,  it  follows  from  equation  (5)  that 
r o  should  also  be  a  constant,  i.e.  center  of  mass  of  the  system  must  remain  fixed 
in  space.  Equation  (5)  provides  the  linear  momentum  conservation  conditions  in  the 
configuration  space. 

Using  Figure  1,  the  position  vector  of  the  end-effector  is  given  as 

r'o  +  lo  +  5Z  +  1P  =  rp(t)  (6) 

1=1 

where  10  is  the  position  vector  of  joint  1  w.r.t.  the  center  of  mass  of  the  manipulator 
base,  1,  is  the  position  vector  of  joint  i  +  1  w.r.t.  joint  i,  lp  is  the  position  vector  of  the 
end-effector  w.r.t.  joint.  n„,  and  rp(t)  is  the  position  vector  of  the  end-effector  from 
the  origin  of  the  inertiat  frame.  The  orientation  of  the  end-effector  may  be  written  in 
a  similar  manner. 

Equations  presented  so  far  are  in  the  vector  form  and  therefore  their  components 

r 

may  be  written  in  any  coordinate  system.  Let  q  =  <71 ,  •  •  • ,  <?n<J  be  a  vec¬ 

tor  of  generalized  coordinates  defining  the  configuration  of  the  system.  Selection  of 
generalised  coordinates  is  not  the  major  issue  in  this  paper.  However,  it  is  worth 
mentioning  that  these  coordinates  must  be  chosen  very  carefully;  otherwise,  it  may 
lead  to  mathematical  singularity,  the  resulting  code  may  be  computationally  ineffi¬ 
cient,  and  in  some  cases  it  may  fail  to  give  the  desired  solution.  In  literature,  one  of 


7 


the  favorite  choices  for  these  coordinates  has  been  to  use  a  set  of  translational  and 
rotational  coordinates  of  the  manipulator  base,  relative  rotation  (or  translation)  of  an 
arm  with  respect  to  its  predecessor  body,  and  orientation  of  the  inertia  wheels  with 
respect  to  the  base.  Base  rotational  coordinates  may  be  represented  by  Euler  angles, 
Euler  like  angles,  a  set  of  three  appropriate  parameters,  or  the  Euler  parameters.  If 
the  Euler  parameters  are  used,  then  one  must  impose  an  additional  Euler  parameter 
normalization  condition. 

In  terms  of  q,  equations  (5)  and  (6)  may  be  given,  respectively,  as 

tfi(q)  =  0  (7) 

<Mq,0  =  0  (8) 

In  equation  (8),  time  t  appears  explicitly  because  of  the  presence  of  rp(f)  in  equation 
(6).  Note  that  if  the  system  is  not  at  a  mathematical  singular  point,  then  the  linear 
and  angular  velocity  terms  are  linearly  related  to  the  time  derivate  of  q.  Based  on 
this  fact,  equations  (3)  and  (4)  and  the  time  derivative  of  equation  (8)  may  be  written 
as 


Jj(q)q  =  0  (9) 

J2(q)q  =  0  (10) 

J3(q,0q  =  X3  (11) 

where  J3(q,  t)  =  (d<f>3j/(d q)  and  X3  =  — (<903)/(<9f).  Matrices  Ji(q),  -Mq),  and 
J3(q,  t)  are  called  the  Jacobian  matrices  associated  with  linear  momentum  (Eq.  3), 
angular  momentum  (Eq.  4),  and  end-effector  kinematic  conditions  (Eq.  6),  :^spec- 
tively.  Equations  (9)  to  (11)  may  be  written  in  a  combined  form  as 

J(q,<)q  =  X  (12) 
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where 


J(q,0  = 


Ji(q) 

J2(q) 

Js(q,0 


is  the  combined  Jacobian  matrix  of  the  system  of  dimension  mj  x  raj,  and 


(13) 


X  = 


(14) 


For  a  redundant  system,  matrix  J  is  in  general  a  rectangular  with  >  mj.  The 
difference  (n^  —  mj),  provides  the  degree  of  redundancy  at  the  velocity  level.  Using 
the  pseudoinverse  of  the  Jacobian  matrix,  the  general  solution  of  equation  (12)  may 
be  written  as  [12] 


q  =  J#X  +  (I- J#J)u  (15) 

where  J *  is  the  pseudoinverse  of  the  Jacobian  matrix  and  u  is  an  arbitrary  unknown 
vector.  The  first  part  of  this  equation  is  the  particular  solution  of  equation  (12).  This 
is  also  the  least  square  solution  of  the  same  equation.  The  second  part  of  equation 
(15)  is  the  complementary  solution  or  the  solution  of  the  homogeneous  equations  [12]. 

From  equations  (12[  and  (15),  the  following  two  observations  may  be  made:  first, 
by  setting  u  =  0  in  ecuation  (15)  one  obtains  only  a  local  minimum  for  q  which 
may  or  may  not  give  a  global  minimum  [13].  Furthermore,  if  ||  X  |j  is  zero,  as  is 
the  case  in  problem  two,  then  equation  (15)  would  lead  to  a  trivial  and  undesirable 
solution.  Second,  if  X  ^  0,  as  may  be  the  case  in  the  first  problem,  then  the  solution 
of  equation  (12)  will  not  lie  in  the  null  space  of  the  Jacobian  matrix  and  the  second 
part  of  equation  (15)  only  will  be  insufficient  to  define  the  solution  of  equation  (12). 
The  present  formulatio  i  directly  works  with  equation  (12)  and  therefore,  it  accounts 
for  both  the  particular  and  the  homogeneous  solutions.  Note  that  by  eliminating 
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equation  (11)  from  equation  (13),  one  obtains  a  set  of  constraints  for  problem  2.  For 
this  reason,  the  second  problem  may  be  considered  as  a  subset  of  the  first  problem. 

In  the  next  section,  we  derive  the  necessary  conditions  for  global  optimum  path 
planning  for  a  space  robotic  manipulator  using  a  variational  approach  [11]. 

3.  Global  Optimization  Formulation 

A  global  optimization  problem  may  be  stated  as  follows: 

rtf 

Min;mize  // =  /  /(q,q,f)df  (16) 

Jto 

subjected  to  the  conditions  in  equation  (12).  Here  If  is  the  objective  functional,  t0 
and  tf  are  the  initial  and  the  final  times,  and  /  is  some  given  function.  Terminal 
times  t0  and  t f  may  be  fixed  or  free.  Some  of  the  choices  for  function  /  considered  in 
the  past  are: 


First  /(q,  q,  0  =  ^qTq,  (17) 

Second  /(q,  q,  f)  =  ^qrA(q,  <)q  (18) 

and 


Third  =  l.  (19) 

/  =  (l/2)qrq  represents  minimization  of  half  of  the  Euclidean  velocity  norm  or  the 
unit-mass-kinetic  energy,  /  =  (l/2)qTAq  represents  minimization  of  the  weighted 
Euclidean  velocity  norm  or  the  kinetic  energy  if  A  is  the  mass  matrix,  and  /  =  1 
represents  minimization  of  time.  In  many  space  applications,  /  may  also  represent 
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objectives  of  obstacle  a/oidance,  singularity  avoidance,  or  minimization  of  the  base 
reaction  force. 

Here  we  consider  that  to  and  tj  are  fixed  and  function  /  is  given  by  equation  (18). 
Using  the  calculus  of  variations,  the  variation  of  equation  (16)  at  the  optimum  point 
may  be  written  as 


61  j  = 


6q(t)dt  + 


'  Of 

dq  (*/) 


6q{tf)  - 


'  df  ~ 
dq(t0) 


8q(t0)  =  0 


(20) 


where  £q(2)  is  the  virtu  al  displacement  of  q(t).  Note  that  the  coefficient  of  <5q(t)  in 
equation  (20)  in  the  present  form  can  not  be  set  equal  to  zero.  This  is  because  all 
components  of  q  are  net  independent.  For  virtual  displacement,  <5q(<)  equation  (12) 
leads  to 


J(q,£)<5q=0  (21) 

Multiplying  equation  (21)  with  AT,  where  A  is  the  vector  of  Lagrange  multipliers,  and 
integrating  the  result  with  respect  to  time  from  t0  to  </,  we  obtain 


ATJ(q,  t)6qdt 


=  0. 


Combining  equations  (20)  and  (22)  we  have 


(22) 


)-  ArJ(q,0 


6q(t)dt  + 


df 

dq{t,) 


6q{tf)  - 


Of 

dq{t0) 


8q(to)  =  0 

(23) 


Using  equation  (23)  and  following  the  discussion  presented  in  [11],  we  obtain  the 
following  Euler- Lagrange  differential  equation 


dq‘ 


&T-769i)T-^t)Tx=a 


and  the  following  terminal  conditions 


(24) 
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— 6q(tb)  =0  tb  =  tQ  or  t}.  (25) 

aqUi). 

Solution  of  equation  (24)  which  satisfies  equations  (7),  (8),  and  (12)  along  with  a 
set  of  boundary  conditions  which  are  consistent  with  equation  (25)  yield  the  desired 
global  optimum  joint  trajectory.  Note  that  some  of  the  concepts  presented  here  are 
similar  to  the  approach  presented  by  Kazerounian  and  Wang  [13].  Their  formulation, 
however,  does  not  account  for  the  nonholonomic  constraints. 

The  above  formulation  is  general  and  applicable  to  a  large  class  of  function  /.  In 
particular,  for  /  given  by  equation  (18),  equations  (24)  and  (25)  reduce  to 

A(q,t)q  +  A(q,t)q-  ^  q  +  J(q,f)TA  =  0,  (26) 

and 

qT(4)A(q(t0),  O^q(^)  =  0  4  =  (o  or  tf  (27) 

where  it  is  assumed  that  A  is  a  symmetric  matrix.  If  not,  then  equation  (18)  can  be 
modified  such  that  the  resulting  matrix  A  is  symmetric.  A  special  case  is  when  A 
is  the  identity  matrix.  This  is  equivalent  to  considering  equation  (17)  for  /.  In  this 
case,  equations  (26)  and  (27)  reduce  to 

q  +  J(q,07’A  =  0,  (28) 

and 

q7(4)<5q(4)  =  0  tb  =  t0ortj  (29) 

Equation  (26)  (or  equation  (28))  is  a  set  of  rid  second  order  differential  equations. 
This  set  is  equivalent  tc  2 rid  first  order  differential  equations.  Alternatively,  one  can 
use  Pontryagin’s  maximum  principle,  which  maximizes  a  certain  Hamiltonian  func¬ 
tion,  to  obtain  a  similar  set  of  2 rid  first  order  equations.  The  calculus  of  variations 
used  here  results  in  simpler  equations  which  are  suitable  for  physical  interpretations 


as  well  as  numerical  and  symbolic  manipulations  [13]. 


4.  Boundary  Conditions 

Equation  (26)  provides  the  necessary  conditions  for  functional  Ij  in  equation  (16) 
to  be  optimum.  Kazerounian  and  Wang  [13]  present  four  different  sets  of  boundary 
conditions  for  a  similar  optimization  problem.  Although,  the  results  of  reference  [13] 
are  applicable  in  this  formulation  also,  they  must  be  used  with  proper  understanding. 
For  the  sake  of  completeness,  we  include  the  results  of  reference  [13]  and  outline  the 
essential  differences. 

In  order  to  obtain  a  proper  set  of  boundary  conditions,  we  return  to  equation 
(27)  (Take  equation  (29)  if  A  =  I).  Note  that  the  generalized  coordinates  are  not  all 
independent.  This  implies  that  the  virtual  displacements  are  al^o  not  all  independent, 
i.e.  all  components  of  5q(ij),  tb  =  to  or  tf, in  equation  (27)  (or  (29))  cannot  be  varied 
freely  and  they  must  satisfy  equation  (21)  at  the  end  points.  Discussion  of  reference 
[13]  cannot  be  applied  here  in  a  straight  forward  manner  because  of  the  presence  of 
nonholonomic  constraints. 

The  necessary  natural  conditions  may  be  evaluated  as  follows:  First,  evaluate 
equation  (21)  for  t  =  tQ,  second,  multiply  the  result  with  nT,  where  fi  is  the  vector  of 
Lagrange  multipliers  of  dimension  mj  x  1,  third,  subtract  the  result  from  equation 
(27)  for  4  =  and  finally  repeat  the  procedure  for  t  =  tj.  This  leads  to 

[qT(4)A(q(t6),f)  -  .uT(ti,)J(q(f6),4)]  <5q(4)  =  0  tb  =  t0ortf  (30) 

Some  remarks  on  \i  wi  1  be  made  in  the  next  section.  In  the  discussion  to  follow, 
consider  that  no  generalized  coordinate  is  specified  at  the  end  points.  In  this  case  the 
coefficient  of  <5q(i{,) ,  (tb  =  to,  t j)  in  equation  (30),  can  be  set  equal  to  zero  for  a  proper 
vector  /i  and  a  proper  set  of  independent  virtual  displacement  (see  reference  [11]). 
Therefore,  for  this  vector  //  and  the  independent  virtual  displacements,  equation  (30) 
after  some  algebra,  leads  to 
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q (h)  =  A  1{q(tb),tb)JT(q(tb),tb)n(tb)  tb  =  t0  or  tf 

Equation  (12)  should  be  satisfied  at  time  tb  also.  This  implies 


(31) 


J(q(4),4)q(4)  =X(4)  (32) 

From  equations  (31)  and  (32),  we  obtain 

n(h)  =  [J^r1Jr]-1X(<6)  (33) 

Substituting  the  value  of  fi  back  in  equation  (31),  we  obtain  the  generalized  velocity 
vector  at  time  tb  as 

q(tb)  =  v4_1Jt[.M-1Jt]-1X(4)  tb  =  t0  or  tf  (34) 

Equation  (34)  is  referred  to  as  the  "natural”  boundary  conditions.  Taking  A  = 
(1/2)/  in  equation  (34)  we  get  the  same  set  of  natural  boundary  conditions  as  that 
of  reference  [13]. 

We  are  now  in  a  position  to  obtain  the  proper  boundary  conditions.  Traditionally, 
equation  (30)  is  used  to  obtain  the  following  four  sets  of  boundary  conditions: 

Set  1.  Generalized  coordinates  are  free  at  both  ends:  In  this  case  the  optimum 
trajectory  must  satisfy  the  natural  boundary  conditions  given  by  equation  (34) 
at  both  end  points. 

Set  2.  Generalized  coordinates  are  given  at  t  =  to  but  not  at  t  =  tj:  In  this  case 
<5q(t0)  =  0  and  therefore  equation  (30)  for  tb  =  t0  is  automatically  satisfied  and 
the  optimum  trajectory  must  meet  the  natural  conditions  at  t  =  tf  only. 

Set  3.  Generalized  coordinates  are  given  at  t  =  tf  but  not  at  t  =  t0:  This  set  is 
identical  to  set  2.  except  that  the  roles  of  to  and  tf  have  interchanged.  In  this 
case  Sq(tf)  =  0  and  equation  (30)  for  tb  =  tf  is  automatically  fulfilled  and  the 
optimum  trajectory  needs  to  agree  with  the  natural  conditions  at  t  =  t0  only. 
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Set  4.  Generalized  coordinates  are  specified  at  both  ends:  In  this  Sq(to)  =  6q(t j)  =  0 
and  equation  (30)  is  automatically  satisfied  for  both  1 1>  =  t0  and  tb  =  tf. 

Sets  1,  2,  and  4  correspond,  respectively,  to  cases  1,  2,  and  4,  of  reference  [13].  In 
the  above  four  sets,  it  is  assumed  that  at  each  end  either  all  (independent)  generalized 
coordinates  are  prescribed  or  all  of  them  (independent  generalized  coordinates)  are 
free.  This  might  not  be  the  case.  In  reality,  part  of  the  (independent)  generalized 
coordinates  may  be  prescribed  and  part  of  them  may  be  free.  Therefore,  there  are 
many  more  possible  combination  sets  which  will  also  satisfy  the  optimality  criteria. 
More  specifically,  if  a  holonomic  system  has  nj  degrees  of  freedom,  then  there  are  22nd 
possible  number  of  sets  which  will  satisfy  the  optimality  criteria.  For  a  nonholonomic 
system,  as  is  the  case  here,  this  number  depends  on  the  number  of  holonomic  and 
nonholonomic  constraints  in  a  complex  way. 

At  this  stage,  it  is  worth  emphasizing  the  following  point:  If  a  system  has  n 
generalized  coordinates,  then  the  possible  number  of  sets  of  independent  coordinates 
is  nC„d,  where  C  is  a  combination  symbol.  These  sets  correspond  to  only  one  set 
in  the  optimality  criteria  discussed  above  because,  for  a  given  set  of  independent 
coordinates,  the  other  coordinates  are  all  fixed.  Therefore,  a  given  set  of  independent 
coordinates  can  be  uniquely  mapped  to  another  set  of  independent  coordinates. 

Among  all  possible  sets,  set  1  stated  above  gives  the  minimum  value  for  the  object 
function,  because  this  sat  imposes  no  restriction  on  the  system.  In  practice,  it  is  also 
possible  to  have  a  set  of  boundary  conditions  which  do  not  satisfy  the  optimality 
criteria.  One  such  set  is  when  the  independent  generalized  coordinates  and  velocities 
are  specified  at  time  t  —  t0.  This  corresponds  to  case  3  of  reference  [13].  This 
set  results  in  an  initial  value  problem.  Since,  in  this  case,  the  natural  boundary 
conditions  at  t  =  tj  are  ignored,  the  resulting  path  is  a  ’’weak”  minimum.  A  strategy 
for  obtaining  strong  minimum  for  this  case  is  given  in  [13]. 

Note  that  the  equations  satisfying  the  optimality  criteria  lead  to  split  boundary 
conditions,  i.e.  two  point  boundary  value  problems.  A  close  form  solution  of  these 
equations  is  generally  not  possible.  Furthermore,  these  equations  can  not  be  solved 
in  a  straight  forward  manner.  A  common  approach  to  this  problem  is  as  follows:  1) 
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estimate  the  initial  conditions,  2)  numerically  integrate  the  differential  and  algebraic 
equations,  3)  use  the  numerical  results  at  t  =  t/  to  update  the  initial  conditions,  and 
4)  repeat  the  steps  2  and  3  until  the  terminal  conditions  are  satisfied  or  the  number 
of  iteration  exceeds  a  prespecified  value.  A  numerical  scheme  to  integrate  the  differ¬ 
ential  and  algebraic  equations  appearing  in  this  formulation  is  given  next. 

5.  Numerical  Integration  of  Differential  and  Algebraic  Equations 

Global  optimum  path  planning  formulation  presented  above  leads  to  a  system  of 
Differential  and  Algebraic  Equations  (DAEs).  Although,  the  formulation  leads  to 
a  two  point  boundary  value  problem,  in  this  section  we  begin  with  some  assumed 
initial  conditions,  i.e.  we  take  some  suitable  values  for  q(t0),  and  q(/0)-  A  method 
to  improve  the  initial  conditions  so  that  the  resulting  solution  satisfies  the  boundary 
conditions  at  both  ends  is  given  in  the  next  section.  The  objective  of  this  section  is 
to  present  a  numerical  scheme  to  advance  the  solution  of  the  DAEs  from  a  time  grid 
point  <,  to  the  next  time  grid  point  tl+\.  For  a  systematic  development  and  ease  of 
reference,  the  DAEs  are  rewritten  below: 

Differential  Equation  (Eq.  (26)): 

A(q,  <)q  +  J(q,  t)T\  =  F( q,  q,  t)  (35) 

Holonomic  Constraints: 

<^i(q) 

4>(q,t)=  =0  (36) 

_  <fo(q,0  _ 

Non-holonomic  Constraints  along  with  the  time  derivative  of  equation  (36): 

J(q,*)q  =  X  (37) 

Vector  F(q,  q,  t)  in  equation  (35)  is  given  as 

•  i  T 

F(q,q,0=J  q-4(q,<)q  (38) 
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Note  that  the  Lagrange  multipliers  A’s  are  unknown.  These  multipliers  may  be  elim¬ 
inated  from  equation  (35)  as  follows:  First,  differentiate  equation  (37)  with  respect 
to  time.  This  leads  to 


J(q,0q  =  x  -  j(q,<)q  (39) 

Second,  substitute  the  value  of  q  from  equation  (35)  into  equation  (39),  and  solve  for 
A.  This  gives 


A=  [jA~lJT]  1(JA~lF-X  +  jq)  (40) 

Finally,  substitute  the  expression  for  A  back  into  equation  (35).  After  some  algebra, 
this  leads  to 


q  =  A~l  F  —  A~lJ  JA~XJT  ( JA~lF  -  X  +  Jq) 


-l 


(41) 


Before  we  proceed  further,  note  that  A  (equation  (40))  at  t  =  t0  (or  tj)  is  different 
from  fj.(t0)  (or  n{tj)).  In  reality,  A  and  /z  are  entirely  two  different  vectors  and  they 
should  be  treated  so. 

The  right  hand  side  of  equation  (41)  contains  no  unknown  terms.  Therefore,  given 
q(<,),  q(£,j,  and  equation  (41)  may  be  used  to  compute  q(f,).  Another  approach 
to  compute  q (<{)  is  as  follows:  Equations  (35)  and  (39)  may  be  written  in  combined 
form  as 


A  iT 

q 

F 

J  0 

A 

X- jq 

(42) 


In  this  equation,  the  square  matrix  and  the  right  hand  vector  can  be  computed 
numerically.  Vectors  c  and  A  then  can  be  solved  using  a  numerical  scheme  such 
as  the  Gaussian  elimination  technique.  Thus,  vector  q  may  be  obtained  using  either 
equation  (41 )  or  equation  (42).  Equation  (42)  has  two  major  advantages  over  equation 
(41).  First,  equation  '42)  preserves  the  sparsity  of  the  matrices  and  therefore  a 
sparse  matrix  code  may  be  used  to  store  the  matrix  elements  and  solve  the  resulting 
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equations,  and  second,  if  the  eigenvalues  of  J  vary  widely  then  equation  (41)  results 
in  ill  condition  matrices  which  cause  numerical  instability.  For  the  same  Jacobian 
matrix,  equation  (42)  leads  to  a  relatively  well  behaved  problem. 

Let  the  solution  for  vector  q  be  written  symbolically  as 

q  =  <?(q,40  (43) 

The  vector  function  g  is  never  computed  explicitly  but  only  numerically.  Equation 
(43)  may  be  reduced  to  a  set  of  first  order  equations  as 


u 

— 

g(  q,u,0 

=  gi(y,t) 

(44) 

q 

u 

Equation  (44)  may  now  be  integrated  using  a  direct  integration  scheme.  This  scheme 
may  lead  to  numerical  problems  because  the  components  of  q  are  not  all  independent 
due  to  the  presence  of  nonlinear  holonomic  and  nonholonomic  constraints.  These 
constraints  nonlinearly  transform  the  small  numerical  integration  errors  which  act 
as  feedback  to  the  system  causing  large  constraint  violation  and  numerical  instabil¬ 
ity.  Furthermore,  since  the  dimension  of  q  may  be  much  larger  than  the  number  of 
independent  coordinates,  these  schemes  require  integration  of  a  much  larger  set  of 
differential  quations. 

The  DAEs  in  this  formulation  are  similar  to  those  in  multi-body  dynamics.  In 
recent  years,  several  methods  have  appeared  in  the  area  of  multibody  dynamics  to 
solve  such  DAEs.  Some  of  these  methods  attempt  to  satisfy  the  constraints  explicitly 
and  some  implicitly,  and  some  take  a  hybrid  approach.  A  brief  review  of  these  methods 
appear  in  reference  [14],  From  theory  of  differential  geometry,  it  is  clear  that  the 
solution  of  these  DAEs  lies  on  the  manifolds  defined  by  these  constraints.  It  is  possible 
to  define  coordinate  systems  on  these  manifolds  and  reduce  the  number  of  differential 
equations  to  its  minimi  m  (see  reference  [14]  and  the  references  therein).  In  this  case, 
the  resulting  generalized  coordinates  are  no  longer  the  original  set  but  the  combination 
of  original  coordinates  with  time  varying  coefficients.  In  this  paper  we  shall  L  esent 
the  numerical  to  solve  the  DAEs  in  terms  of  a  set  of  coordinates  which  is  a  subset  of 
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the  original  set  of  generalized  coordinates. 

In  order  to  discuss  the  numerical  scheme  to  follow,  let  q/i(<«)  and  q/2(U)  be  the 
vectors  of  independent  generalized  coordinates  and  velocities  at  any  time  £,.  These 
two  vectors  completely  define  the  state  of  the  system,  because  they  may  be  used  to 
solve  the  dependent  generalized  coordinates  and  velocities.  Note  that  the  dimensions 
of  these  two  vectors  are  not  the  same  because  of  the  presence  of  nonholonomic  con¬ 
straints.  Once  these  vectors  have  been  identified,  the  numerical  scheme  to  solve  the 
DAEs  may  be  given  as  follows: 

Numerical  Algorithm: 

Step  1.  Use  q/i(£,)  and  equation  (36)  to  solve  the  unknown  dependent  coordinate. 
This  requires  solution  of  a  set  of  nonlinear  equations.  Newton- Raphson  method 
and  its  variants  may  be  used  for  this  purpose. 

Step  2.  Use  q/2(£,),  equation  (37)  and  the  vector  q(£t)  obtained  in  step  1  to  solve 
for  the  dependent  velocities.  This  requires  solution  of  a  set  of  linear  equations. 
Any  of  the  many  schemes  such  as  Gaussian  elimination,  etc.  may  be  used  for 
this  purpose. 

Step  3.  Solve  for  q(£i)  using  either  equation  (41)  or  (42)  and  identify  vector  q/2(£i) 
from  this  vector. 

Step  4.  Use  vectors  q/2(£i)  and  q/2(£i)  in  an  integration  subroutine  to  obtain  vectors 
q/2(U+i)  and  q/i(£,+i). 

Step  5.  Repeat  steps  .  to  4  until  final  time  has  reached. 

Steps  1  and  2  insure  that  the  kinematic  conditions  are  satisfied.  Also  note  that 
the  number  of  equations  that  needs  to  be  integrated  is  much  less  than  the  number  of 
differential  equations  in  (44).  In  this  approach,  however,  one  must  solve  additional 
linear  and  nonlinear  equations. 
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6.  Iterative  Shooting  Method 


A  numerical  method  to  solve  the  DAEs  was  presented  in  the  previous  section. 
As  shown  above,  the  problem  addressed  here  leads  to  a  two-point  boundary  value 
problem.  The  assumed  initial  guess  in  the  previous  section  generally  does  not  satisfy 
the  boundary  conditions.  The  guess  may  be  improved  using  a  shooting  method.  The 
basic  idea  of  the  shooting  method  is  as  follows:  First,  start  with  an  initial  guess  and 
find  the  response  at  the  final  time,  second,  use  some  numerical  scheme  to  estimate 
the  initial  guess,  and  third,  repeat  the  first  and  the  second  steps  until  the  boundary 
conditions  are  satisfied  or  the  iteration  limit  is  reached. 

In  this  section  we  consider  the  Secant  method  to  improve  the  initial  guess.  This 
is  a  well  established  technique  and  the  detailed  numerical  algorithm  for  this  may 
be  obtained  in  many  of  the  standard  numerical  analysis  literature.  Here,  we  briefly 
outline  the  method.  Let  the  two-point-boundary  value  problem  be  as  follows:  Find 
the  response  of  the  system 


y  =  /(y,y,0  (45) 

which  satisfies  the  following  boundary  conditions:  y(i0)  =  yo  and  y (tj)  =  y/.  In  this 
method  two  initial  guesses  for  y(to)  are  considered.  The  Secant  method  for  improving 
the  initial  guess  for  y(t0)  =  y0  is  given  as  follows: 


yk(to)  =  yk-'(t0)- 


(yk  \tj)  -ys){yk~x{tQ)  -yk~\t0)) 


(46) 


yk~l(h )  -  yk~2{tf) 

where  k  is  the  iteration  number.  This  scheme  generates  a  sequence  of  yk(to)  such 
that 


lim  yk{tf)  =  y7 

k—+oo 

Thus,  the  initial  guess  that  satisfies  the  two-point-boundary  conditions  may  be 
obtained. 
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7.  Numerical  Results 

In  order  to  demonstrate  the  feasibility  of  the  formulation,  we  consider  a  two- 
dimensional  space  manipulator  consisting  of  a  base  (body  0)  and  three  arms  (bodies  1, 
2,  and  3)  as  shown  in  Figure  2.  The  link  lengths  are  Li  =  Z2  =  7.0m,  and  —  4.0m. 
Initially,  the  inertia  properties  considered  are  mo  =  100%.,  m\  =  m2  =  7.0%., 
m3  =  4.0%.,  Iq  =  1.0(%  —  m2),  Ii  =  I2  =  28.583 (%  —  m2),  and  =  5.333(%  —  m1 ). 
The  tip  of  the  vector  is  required  to  move  along  a  circular  trajectory  defin'd  as 

-  nr 

xJt)  xc  +  rcos(a ) 

X  =  PW  =  V  ;  (47) 

yP(t)  [  yc  +  T’sm(a) 

where  (xc,yc)  is  the  center  of  the  circle,  r  its  radius,  and  a  a  specified  time  dependent 
parameter.  The  center  is  considered  it  ic  =  11.0m  and  yc  =  6.0m,  and  the  radius 
r  =  3.0m.  a(t)  considered  is 

,  2; r  .  ,27rT  , 

af^  =  y-«n(— )  (48) 

The  objective  is  to  find  au  optima’  trajectory  for  a  set  of  specified  initial  generalized 
r  ordinates.  Also,  the  initial  and  the  final  velocities  should  be  zero.  Furthermore, 
o  robot  ba^e  is  allowed  to  translate,  but  not  to  rotate.  The  initial  conditions,  which 
.re  consistent  with  the  constraints,  are  xo  =  —1.0213,  j/o  =  —0.7124,  0\  =  0.93478, 
02  =  0.0064397,  and  03  =  0.262,  and  all  initial  velocities  are  zero. 

The  numerical  results  of  the  system  are  given  below.  Figure  3  shows  the  dis¬ 
placements  Axo  and  A  jo  as  a  function  of  time.  From  tins  figure,  it  is  clear  that  the 
motion  of  the  arms  considerably  deflects  the  center  of  mass  of  the  base.  The  changes 
in  the  orientations  of  the  arms  as  a  function  of  time  are  shown  in  Figure  4.  Figure 
5  shows  configurations  of  the  space  manipulator  at  various  times.  Numerical  results 
for  this  response  were  also  studied  using  an  animation  program.  It  is  observed  that 
the  end-effector  traces  a  perfect  circle  and  the  response  of  the  system  is  smooth  and 
continuous. 

In  a  further  investigation,  the  mass  and  the  inertia  of  each  arm  were  increased, 


21 


Figure  2.  Configuration  of  a  Two  Dimensional  Space  Manipulator 
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Figure  3.  Time  Response  of  the  Center  of  Mass  of  the  Base 


23 


Figure  4.  Orientations  of  the  Manipulator  Arms  as  a  Function  of  Time 
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Figure  5.  Configuration  of  the  Manipulator  at  Various  Points 
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first  by  three  times  and  second  by  six  times.  Corresponding  changes  in  time  response 
for  Axo?  At/Oi  A$i,  A#2i  and  A03  are  shown  in  Figures  6  to  10.  It  can  be  seen  that 
the  variation  in  the  response  of  the  base  increases  with  increase  in  inertia.  This  is  in 
agreement  with  the  intution. 

8.  Conclusions: 

An  optimum  path  planning  formulation  for  a  free  floating  space  robotic  manipula¬ 
tor  based  on  a  variational  approach  has  been  presented.  The  formulation  shows  that 
the  moment  conservation  conditions  result  in  a  system  of  holonomic  and  nonholo- 
nomic  constraints.  A  method  to  incorporate  these  constraints  has  also  been  given. 
This  leads  to  a  system  of  differential  and  algebraic  equations,  and  two-point  bound¬ 
ary  conditions  if  the  initial  conditions  are  not  specified.  A  numerical  algorithm  to 
solve  the  resulting  DAEs  has  been  proposed.  Also,  a  numerical  scheme  to  solve  the 
two  point  boundary  value  problem  has  been  outlined.  The  formulation  has  been  used 
to  obtain  the  response  of  a  two-dimensional  redundant  space  robot.  The  numeri¬ 
cal  results  show  that  the  motion  of  the  arms  can  considerably  affect  the  response 
of  the  base.  The  approach  is  of  significance  in  many  space  applications  since  most 
current  space  robots  are  redundant.  This  redundancy  provides  the  necessary  dexter¬ 
ity  required  for  extra-vehicular  activity  or  avoidance  of  potential  obstacles  in  space 
stations. 
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Figure  9.  Effect  of  Change  of  Mass  on  theta  ^t) 
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Figure  10.  Effect  of  Change  of  Mass  on  theta  (t) 
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APPENDIX 

Equations  (7)  to  (11)  are  the  main  equations  in  the  derivation  for  global  optimum 
path  planning  for  a  redundant  space  robotic  manipulator.  Once  these  equations  are 
known,  the  other  equations  may  be  obtained  using  the  formulation  presented  in  this 
paper.  In  this  section  we  demonstrate  these  equations  for  a  two  dimensional  space 
manipulator  system. 

In  order  to  accomplish  the  above  stated  objective,  consider  a  two  dimensional 
space  manipulator  consisting  of  one  base  (body  0)  and  three  arms  (bodies  1,  2,  and 
3)  as  shown  in  figure  2.  The  configuration  of  this  system  is  defined  using  two  sets 
of  frames,  namely  an  inertial  frame  and  four  body  frames.  The  origin  of  each  body 
frame  is  rigidly  attached  to  the  center  of  mass  of  the  respective  body.  For  simplicity, 
it  is  assumed  that  the  total  linear  and  angular  momentums  of  the  system  are  zero 
and  the  origin  of  the  inertial  frame  coincides  with  the  center  of  mass  of  the  overall 
system.  For  this  system,  the  vector  of  generalized  coordinates  q  taken  here  is 


q  = 


02 


03  *o  yo 


(A  1) 


where  0,  is  the  orientat'  on  of  the  x-axis  of  body  frame  i  with  respect  to  the  x-axis  of 
the  inertial  frame,  and  (xoij/o)  is  the  location  of  center  of  mass  with  respect  to  the 
inertial  frame.  Let  m,  and  /,  be  the  mass  and  moment  of  inertia  of  body  i  and  m  be 
the  total  mass  of  the  system;  i.e. 


m  =  m0  +  mi  +  m2  +  m3  (.42) 

The  linear  and  angular  momentum  conservation  conditions  for  this  system  lead  to 

3 

]Tm,rt  =  0  (.43) 

1=0 

3 

E(/.0t  +  m'r>  x  r.)  =  0  (.44) 

i=0 

where  r,  is  the  position  vector  of  center  of  mass  of  body  i,  and  the  period  on  (*) 
denoted  total  time  derivate  of  (*).  Vectors  r,  ( i  =  0, 1,2,3)  may  be  written  in  terms 
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of  q  directly  from  figure  2.  Substituting  the  expression  for  r,  in  equations  (A3)  and 

(A4),  we  obtain  equations  (8)  and  (9).  For  this  case,  the  components  of  matrices  Ji 

and  J2  are  as  follows: 

[Ji]n  =  +  m2  +  m3)  *  Li*  sin(0x) 

(A5) 

771 2 

[J  1  ]  1 2  =  —  (~ — V  m?>)  *  L 2  *  sin(92) 

(-46) 

[Jt]i3  =  *  £3  *  sm{63) 

(.47) 

[  J 1  ]  14  =  m 

(.48) 

C-l 

1— 1 

1— • 

CM 

II 

O 

(.49) 

[Ji]ie  =  (m  -  m0)  *  (-£  *  sin(90)  —  V  *  cos(00)) 

(.410) 

[J  1(21  =  (-^ — F  tti2  +  m3)  *  L\  *  cos(Q\ ) 

(.411) 

1712 

[Ji]22  =  (“2"  +  m3)  *  L2*  cos(92) 

(A  12) 

[ J 1  ] 23  =  ~  *  L3  *  COs(03) 

(A13) 

[  J 1  ]  2 1  “  0 

(-414) 

[  J 1  j  25  =  m 

(A  15) 

[ J i ] 26  =  (m  -  m0)  *  (£  *  cos(0o)  —  V  *  sin(90))  (.416) 
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[J2]n  =  h  4-  4-  m2  +  m3)Li(x0cos(9i)  4  y0sin(9i)  +  £cos(9x  -  90)  +  r}sin(Ox  -  90)) 

+(— — h  m2  4  tti3)L \  4  (— — I-  m3)L\L2Cos{9\  —  #2)  4 — —L\L3cos(Qi  —  03)  (417) 


[J 2 J 12  —  I2  4-  (— — f  m3) L,2(xocos(92.)  4  yosin(02)  +  <f cos(#2  —  $0)  4  rjsin(0 2  —  $o))4 
(— — f  m3)LiL2cos(92  —  9i )  4  (— — f-  m3)L\  4 — —  L,2L3cos(92  —  93)  (.418) 


[J2]i3  =  h  +  —  L3((x0cos(93)  4  y0sin(93)  4  {cos(93  -  90)  4-  rjsin(93  -  90) 


+  L\cos{93  —  9i)  +  Li2Cos(93  —  #2)  4  ~^L3) 


[J2]i4  =  -my0  -  (m  -  mo)(£sm(0o)  4  tjcos(90))  -  (—  4  m2  4  m3)Lxsin{Ox) 

mo  m-x 

-(—4 -  m3)L2sin(92) - —  L3sin{93)  (4 


[ J2]  15  =  mx 0  4  (m  -  mo)(£co.s(0o)  -  r]sin(0o))  4  (  —  4  m2  4-  m3)Lxcos(9x  \ 
+(^~  4  m3)L2cos(92)  4  ^ L3cos{93)  ( 


[J2]ie  =  Io+(Tn-mn)(Z2+y2)+((m-mo)xo+(T?Y+rn2+m3)Llcos(0i)+(T^+m3)L2Cos(O2 


+  L3cos{93))(£cos(00)  -  r}sin(90))  4  ((m  -  m0)y0  4  4  m2  4  m3)Lisin(Oi) 4 

4  m3)L2sin(02)  4  ^ L3sin(93))(£sin(90)  4-  ycos{0o))  (.422) 

Equations  (A20)  and  (A21)  follow  from  linear  momentum  conservation  conditions. 
Integrating  equation  (A3)  and  using  the  fact  that  the  center  of  mass  of  the  system 
lies  at  (0,0),  we  obtain  equation  (7)  as 


34 


<Mq)  = 

where  functions  <^u(q)  and  <t> i2(q)  are  given  as 


- 

<Mq) 

r  ■) 
0 

<^i2(q) 

0 

(A23) 


<^n(q)  =  mio  +  (m  —  mo)(^cos(0o)  -  rjsin(90))  +  (~-  +  m2  4-  Tn.3)LiCos(9x)+ 


(~  +  m3)l2cos(02)  +  ^~L3cos(9 3)  =  0 


(.4.24 1 


^nCq)  =  mj/o  +  (m  -  mo)(£sin(0o)  +  t/cos(^0))  +  (y-  +  m2  +  m3)Z,1sfn(01)  + 


(-y  +  m3)I2^n(02)  +  -y L3sin{93 )  =  0 


(A25) 

T 

be 


Equations  (A24)  and  (A25)  lead  to  [J2]i4  =  [J2]u  =  0.  Let  X3  =  Xp(t)  yp{t) 
the  specified  trajectory  of  the  end  effector.  Then  the  trajectory  constraint  (equations 
(8))  is  given  as 


<£3i(q,0 

0 

_  <Mq,0 

0 

4>3(q,0  = 


where  functions  03i(q,  i)  and  <^32(q,  t)  are  given  by  the  following  equations: 


(.426) 


<t>3  i(q,0  =  xo+(,*cos{9G)-Ti*siTi{0o)  +  Li*cos{9l)  +  L2*cos{92)  +  L3*cos{93)-xp(t)  =  0 

(.427) 


<^32(q,  t)  =  Vo  +  (,  *  sin(Oo)  +  TJ*  cos(Qq)  +  sin{9\  )  +  L2*  sin(92)  +  L3*  sin(03)  -  yp(t) 


In  this  study,  the  trajectory  is  a  circle  specified  as 


(.428) 


xp(/.) 

_  uP{t) 

~ 

xc  +  rcos(a{t ))  yc  +  rsin(a(t)) 

(.429) 
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where  (xc,  yc)  and  r  are,  respectively,  the  coordinate  of  the  center  and  the  radius  of 
the  circle,  and  a(t)  is  a  specified  time  dependent  parameter.  In  this  study,  a(t)  is 
considered  as 


.  .  2ir  .  2 nt 

=  y  ~ 

(.430) 

Time  derivative  of  equations  (A27)  and  (A28)  gives  equation  (11). 

of  the  Jacobian  matrix  J3  are  as  follows: 

The  components 

[ J 3)11  =  ~L\  *  sin(di) 

(-431) 

[J3]  1 2  =  —  L?  *  sin(6 2) 

(.432) 

[Jslia  =  -L3  *  sin(d3) 

(.433) 

[J3]  14  =  1 

(-434) 

[*^3]  15  =  0 

(A35) 

[J 3]  16  =  -£  *  sin(Oo)  -  rj  *  cos(0o ) 

(.436) 

[ J3] 2 1  —  L\  *  cos(Ox) 

(.437) 

[J  3]  22  =  Li*  cps($i) 

(.438) 

[J 3] 23  =  L3  *  COS(03 ) 

(.439) 

[  J3]  24  =  0 

(.440) 
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[•J  3]2S  =  1 


M41) 


[^26  =  £  *  cos(90)  -  t)  v  sin(9o) 
and  the  vector  X3  is  given  as 

—rsin(a) 
rcos(a) 

Thus,  equations  (7)  to  (9)  are  known  for  the  system. 


•  2tt  2tt< 

X3  =  t(1_cos(~ 


(A42) 


(.443) 
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